Assignment 1
Implement a Basic Image Processing Pipeline
Honi Selahaddin (20226089) | honis@kaist.ac.kr
Initials
The original image, im0, in 2856x4290 resolution is imported to my workspace in uint16 type. We cast this data type to double which increase 'Bytes' value from 24504480 to 98017920.
im0 = imread('data/banana_slug.tiff');
whos im0
im0 = double(im0);
whos im0
Linearization
This routine shifts all pixel values left by val_min. Then, values below zero are moved to zero, above val_max are moved to val_max as a clipping operation which makes intensity range [0, val_max]. Finally, a normalization is applied to map values into [0, 1]. You can see the histogram of linearized image and visualization via increased brightness below.
val_min = 2047;
val_max = 15000;
im = im0 - val_min;
im(im<0) = 0;
im(im>val_max) = val_max;
im = im/val_max;
figure
subplot(2,1,1)
imhist(im)
title('Histogram of Linearized Image')
subplot(2,1,2)
imshow(min(1, im * 5));
title('The linear RAW image (brightness increased by 5)')
Identifying the Correct Bayer Pattern
In this section, my idea to identify the correct Bayer pattern is cropping 2x2 pixel areas from color-known (or predictable) regions in the image. Top-left part of the image frame keeps green grass (contextually) which requires green buckets to have higher intensity. Below figure shows that intensities of 2nd and 3rd pixels of top-left 2x2 is much more higher than the others. Thus, these pixels are Green. Moreover, bottom-left of the image has a blue-ish color (woman's striped shirt) which leds 4th pixel of the Bayer pattern to have more intensity. At the end, we can conclude that this Bayer pattern is RGGB as specified for Canon EOS T3 Rebel camera.
roi = 2; brightness_param = 40;
figure;
subplot(1,2,1)
imshow(min(1, im(1:roi, 1:roi) * brightness_param));
title('2x2 RoI (top-left) (must be green)')
subplot(1,2,2)
imshow(min(1, im(end-roi+1:end, 1:roi) * brightness_param));
title('2x2 RoI (bottom-left) (must be blue)')
We separate RGB channels from Bayer patterns by simple indexing followed by concatenating operation.
byr1 = im(1:2:end, 1:2:end); %R
byr2 = im(1:2:end, 2:2:end); %G
byr3 = im(2:2:end, 1:2:end); %G
byr4 = im(2:2:end, 2:2:end); %B
im_r = byr1;
im_g = byr3;
im_b = byr4;
im_rgb = cat(3, im_r, im_g, im_b);
figure
imshow(min(1, im_rgb * 5));
title('RGB Construction via Bayer Patterns (brightness increased by 5)')
White Balancing
Human visual system has chromatic adaptation to perceive white colors correctly in varying luminance environments; yet, there is no camera perception. Therefore, we need to balance the white instensity. There are two main assumptions for white balancing: grey-world and white-world.
In grey-world assumption, we normalize each channel by its channel-wise average. Then, normalize all channels by green average. It is good to remember that green channel is more important to understand real intensity for human-eye. In white-world assumption, on the other hand, we use channel-maximums instead of averages.
% Grey-world assumption
im_r_avg = mean(im_r);
im_g_avg = mean(im_g);
im_b_avg = mean(im_b);
im_r_wb_g = im_r * (im_g_avg/im_r_avg);
im_g_wb_g = im_g;
im_b_wb_g = im_b * (im_g_avg/im_b_avg);
% White-world assumption
im_r_max = max(max(im_r));
im_g_max = max(max(im_g));
im_b_max = max(max(im_b));
im_r_wb_w = im_r * (im_g_max/im_r_max);
im_g_wb_w = im_g;
im_b_wb_w = im_b * (im_g_max/im_b_max);
After each channel is white-balanced, we concat these channels to obtain RGB image. You can see the difference of grey-world and white-world assumption results in below figure. (in both frames brightness increased by 5) The colors in grey-world assumption looks more natural; thus, this progress is picked for further pipeline.
im_rgb_wb_g = cat(3, im_r_wb_g, im_g_wb_g, im_b_wb_g); %grey
im_rgb_wb_w = cat(3, im_r_wb_w, im_g_wb_w, im_b_wb_w); %white
figure
subplot(2,1,1)
imshow(min(1, im_rgb_wb_g * 5));
title('WB image under grey-world assumption')
subplot(2,1,2)
imshow(min(1, im_rgb_wb_w * 5));
title('WB image under white-world assumption')
Demosaicing
In this step, we simply apply built-in interp2() function for each channel for demosaicing.
im_r_wb_g_interp = interp2(im_r_wb_g);
im_g_wb_g_interp = interp2(im_g_wb_g);
im_b_wb_g_interp = interp2(im_b_wb_g);
im_rgb_wb_g_interp = cat(3, im_r_wb_g_interp, im_g_wb_g_interp, im_b_wb_g_interp);
figure;
imshow(min(1, im_rgb_wb_g_interp * 5));
title('Bilinear Interpolated Image (brightness increased by 5)')
Gamma Correction
After the shifting and scaling operation in linearization, most pixel intensities got closer to zero and the image darkened.
Now, we convert white-balanced and demosaiced image into gray-scale via rgb2gray() function. This function converts RGB values to grayscale values by forming a weighted sum of the R, G, and B components:
0.2989 * R + 0.5870 * G + 0.1140 * B
Based on gray-scale intensity we perform gamma-correction according to sRGB Standard.
f(u) = c ⋅ u, 0 ≤ u < d
f(u) = a ⋅ uɣ + b, u ≥ d,
where u represents a color value with these parameters:
a = 1.055
b = –0.055
c = 12.92
d = 0.0031308
ɣ = 1/2.4
(from MATLAB Documentation)
im_gray = rgb2gray(im_rgb_wb_g_interp);
if im_gray <= 0.0031308
im_rgb_wb_g_interp_gamm = im_rgb_wb_g_interp * 12.92;
else
im_rgb_wb_g_interp_gamm = im_rgb_wb_g_interp.^(1/2.4) * (1+0.055) - 0.055;
end
figure;
imshow(im_rgb_wb_g_interp_gamm);
title('Brightness Adjusted and Gamma-Corrected Image')
Compression
Image processing pipeline is now complete. The image size without compression is 13.2MB in PNG format. It is time to compress the image for storage efficiency. We can get 2.41MB size with 95% quality in JPEG which has 5.5 compression ratio. Whereas 50% quality in JPEG provides 23.6 compression ratio. You can find other experimental results in following figure. Image starts to lose contextual information after compression ratio of 50.
In common practice people pick the output of compression ratio 10 but my preference in this application is 50% quality JPEG due to its small size 0.56MB while preserving the quality.
imwrite(im_rgb_wb_g_interp_gamm, 'data/im_out.png');
imwrite(im_rgb_wb_g_interp_gamm, 'data/im_out_95.jpeg', 'quality', 95);
imwrite(im_rgb_wb_g_interp_gamm, 'data/im_out_50.jpeg', 'quality', 50);
imwrite(im_rgb_wb_g_interp_gamm, 'data/im_out_10.jpeg', 'quality', 10);
imwrite(im_rgb_wb_g_interp_gamm, 'data/im_out_05.jpeg', 'quality', 5);
imwrite(im_rgb_wb_g_interp_gamm, 'data/im_out_01.jpeg', 'quality', 1);
impng = imread('data/im_out.png');
im95 = imread('data/im_out_95.jpeg');
im50 = imread('data/im_out_50.jpeg');
im10 = imread('data/im_out_10.jpeg');
im05 = imread('data/im_out_05.jpeg');
im01 = imread('data/im_out_01.jpeg');
figure;
subplot(3,2,1)
imshow(impng)
title('PNG 13.2MB')
subplot(3,2,2)
imshow(im95)
title('JPEG95 2.41MB | 5.5')
subplot(3,2,3)
imshow(im50)
title('JPEG50 0.56MB | 23.6')
subplot(3,2,4)
imshow(im10)
title('JPEG10 0.26MB | 50.8')
subplot(3,2,5)
imshow(im05)
title('JPEG5 0.21MB | 62.8')
subplot(3,2,6)
imshow(im01)
title('JPEG1 0.19MB | 69.5')
